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Abstract 

We have studied the relativistic beamed outflow proposed to occur in the 

^ ■ collapsar model of gamma-ray bursts. A jet forms as a consequence of an 

\Q ' assumed energy deposition of ~ 10'''" — 10''^ erg/s within a 30° cone around 

\^ , the rotation axis of the progenitor star. The generated jet flow is strongly 

^^ ' beamed (^ few degrees) and reaches the surface of the stellar progenitor 

^^ ' (r Ri 3 10^" cm) intact. At break-out the maximum Lorentz factor of the 

^"^ , jet flow is about 33. Simulations have been performed with the GENESIS 

^1^ ' multi-dimensional relativistic hydrodynamic code. 

CLi' 1 Motivation and numerical setup 

^ ■ Various catastrophic collapse events have been proposed to explain the energies 

"^ I released in a gamma-ray burst (GRB) including compact binary system mergers 

C^ ■ [H |l^, collapsars [18| and hypernovae [Q. These models all rely on a common 

engine, namely a stellar mass black hole (BH) which accretes several solar masses 
of matter from a disk (formed during a merger or by a non-spherical collapse). 
A fraction of the gravitational binding energy released by accretion is converted 
j^ ■ into a pair fireball. Provided the baryon load of the fireball is not too large, the 

baryons are accelerated together with the e""" e~ pairs to ultra-relativistic speeds 
(Lorentz factors > 10^; H). The existence of such relativistic flows is supported 
by radio observations of GRB 980425 @. 

The dynamics of spherically symmetric relativistic fireballs has been studied by 
several authors by means of ID Lagrangian hydrodynamic simulations (e.g., |12[). 
It has been argued that the rapid temporal decay of several GRB afterglows is 
more consistent with the evolution of a relativistic jet after it slows down and 
spreads laterally than with a spherical blast wave P|. The lack of a significant 
radio afterglow in GRB 990123 provides independent evidence for jet-like geometry 
[ pO| . Motivated by these observations and by the collapsar model of ||ll| , we have 
simulated the propagation of jets from collapsars using relativistic hydrodynamics. 
In pl| the continued evolution of rotating helium stars, whose iron core collapse 
does not produce a successful outgoing shock but instead forms a BH surrounded 
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by a compact accretion disk, has been explored. Assuming that the efficiency of 
energy deposition by i^i?-annihilation or, e.g.,magneto-hydrodynamic processes is 
higher in the polar regions, [O obtained relativistic jets along the rotation axis, 
which remained highly focused, and capable of penetrating the star. However, as 
these simulations were performed with a Newtonian hydrodynamic code, appre- 
ciably superluminal speeds in the jet flow were obtained. 

We have performed axisymmetric relativistic simulations of jets from coUap- 
sars starting from Model 14A of [^. The simulations have been performed with 
GENESIS a multidimensional relativistic hydrodynamic code (based on Godunov- 
type schemes) developed by Q] using 2D spherical coordinates (r, 9) . GENESIS 
employs a 3th order explicit Runge-Kutta method |lq| to advance in the time 
the relativistic Euler equations written in conservation form. High spatial order is 
provided by a PPM reconstruction Q that sets up the values of the physical vari- 
ables in order to solve linearized Riemann problems at every cell interface (using 
the Marquina's flux formula y|). 

The innermost 2.03 Mq representing the iron core were removed from the he- 
lium star model by introducing an inner boundary at a radius of 200 km. When the 
central BH has acquired a mass of 3.762 Mq, we mapped the model to our com- 
putational grid. In the i — direction the computational grid consists of 200 zones 
spaced logarithmically between the inner boundary and the surface of the helium 
star at _R» = 2.98 x 10^" cm. Assuming equatorial-plane symmetry we use four 
different zonings in the angular direction: 44, 90 and 180 uniform zones (i.e., 2°, 1° 
and 0.5° angular resolution), and 100 nonuniform zones covering the polar region 
0° < 6* < 30° with 60 equidistant zones (0.5° resolution) and the remaining 40 
zones being logarithmically distributed between 30° < 9 < 90°. 

The gravitational field of the BH is described by the static Schwarzschild met- 
ric, neglecting the effects due to self-gravity of the star. We used the EOS of |I^ 
which includes the contribution of non-relativistic nucleons treated as a mixture 
of Boltzmann gases, and radiation, as well as an approximate correction due to 
pairs e+e^. Full ionization and non-degeneracy of the electrons is assumed. We 
advect (i.e., we do not solve additional Riemann problems for each component) 
nine non-reacting nuclear species which are present in the initial model. 

In a consistent collapsar model the jet will be launched by any physical process 
which gives rise to a local deposition of energy and/or momentum. We mimic 
this process by depositing energy at a constant rate, E, within a 30° cone around 
the rotation axis of the progenitor star. In radial direction the deposition region 
extends from the inner boundary to a radius of 6 x 10 cm. We consider two cases 
that bracket the expected E of the collapsar models: 10^" erg/s, and 10^^ erg/s. 

2 Results 

Low energy deposition rate (Model A). Using a constant E = 10^° erg/s a 
relativistic jet forms within a fraction of a second and starts to propagate along 
the rotation axis (Fig. 1). The jet exhibits all the typical morphological elements 
[0: a terminal bow shock, a narrow cocoon, a contact discontinuity separating ste- 




Figure 1: Coloured contour maps of the logarithm of the rest-mass density (six top 
panels) and the Lorentz factor for model A at different evolution times. Note the change 
in the scale between left and right panels. 



liar and jet matter, and a hot spot. The propagation of the jet is unsteady, 
because of density inhomogeneities in the star. The Lorentz factor of the jet, W, 
increases non-monotonically with time, while the density drops to ~ 10~^ gr/cm^. 
The density profile shows large variations (up to a factor of 100) due to internal 
shocks. The mean density in the jet is ~ 10~^ — 1 g/cm'^. 

Some of the internal shocks are biconical and recoUimate the beam. These 
shocks develop during the jet's propagation and may provide the "internal shocks" 
proposed to explain the observed gamma-ray emission [[7|. A particularly strong 
recollimation shock forms during the early stages of the evolution, followed by a 
strong rarefaction that causes the largest acceleration of the beam material giving 
rise to a maximum in W. When the jet encounters a region along the axis where 
the density gradient is positive the jet's head is decelerated, while the a central 
channel in the beam is cleaned by outflow into the cocoon through the head. This 
leads to an acceleration of the beam. The combination of both effects (deceleration 
of the head and beam acceleration) increases the strength of the internal shocks. 

The relativistic treatment of the hydrodynamics leads to an overall qualita- 
tively similar evolution than in pT[ (formation of a jet), being, however, a quan- 
titatively very different. We find that the results strongly depend on the angular 
resolution, and the minimum acceptable one is 0.5° (at least near the axis). At this 
resolution we find W^ax '^15 — 20 (at shock break-out) at a radius ~ 8 x 10^ cm. 
Within the uncertainties of the jet mass determinations due to finite zoning and 
the lack of a precise numerical criterium to identify jet matter, the baryon load, 
77, seems to decrease with increasing resolution. In the highest resolution run we 
find 77 ~ 1.3 ± 1.2 at shock break-out (see also Sect. 4). 

High energy deposition rate (Model B). Enhancing E hy a tenfold {E = 
10^-^erg/s), the jet flow reaches larger values of Wmax- We observe transients 
during which Wmax becomes as large as 40 (Wmax = 33.3 at shock breakout). The 
jet propagates faster than in model A. The time required to reach the surface of 
the star is 2.27s instead of 3.35 s. The opening angle of the jet at shock breakout 
is ~ 10°, i.e., the jets is less collimated than model A. The strong recollimation 
shock present in the model A is not so evident here. Instead, several biconical 
shocks are observed, and W near the head of the jet is larger (~ 22 in the final 
model) because, due to the larger E, the central funnel is evacuated faster, and 
because the mean density in the jet is 5 times smaller than in model A [rj being 
twice as large). 

Evolution after shock breakout. After reaching the stellar surface the 
relativistic jet propagates through a medium of decreasing density continuously 
releasing energy into a medium whose pressure is negligible compared to that in 
the jet cavity, and whose density is (initially) of the same order as that of the 
jet. These are jump conditions that generate a strong blast wave. The external 
density gradient determines whether the shock will accelerate or decelerate with 
time (Pql). In order to satisfy the conditions for accelerating shocks (p5|), we have 
generated a Gaussian atmosphere matching an external uniform medium. We use 
models A and B to simulate the evolution after shock breakout. The computational 
domain is extended for this purpose to a radius of Rt — 7.6 x lO^*' cm. The jet 
reaches Rt (from the stellar surface) after 1.8s in both models, i.e., the mean 
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Figure 2: Evolution of the axial and lateral sizes of the jet cavity during the post- 
breakout epoch. Time is measured with respect to the breakout time for each model. 



propagation velocity is ^ 0.85c (almost three times larger than that inside the 
star) . 

The evolution after shock breakout can be distinguished into three epochs (see 
Figs.H and 0), which are related with (i) the external thermodynamical gradients 
and (ii) the importance of the axial momentum flux relative to the pressure into the 
jet cavity. Both effects determine the shape of the expanding bubble -prolate- (see 
Figs.j^ and g) during the post-breakout evolution. However, when the jet reaches 
the uniform part of the circumstellar environment, the shape changes appreciably, 
because the sideways expansion is faster. We have not followed the evolution long 
enough to see what happens when most of the bubble has reached the uniform part 
of the environment. Nevertheless, we can infer from Fig. 2 than the widening rate 
reduces with time in a way similar to what has happened to the axial expansion. 
At latter times most of the bubble is inside the uniform medium, and the bubble 
will eventually be pressure driven. Hence a isotropic expansion is expected. 

After shock breakout there are transients in which Wmax becomes almost 50 in 
some parts of the beam, Wmax is again obtained behind the strongest recoUimation 
shock. The Lorentz factor near the boundary of the cavity blown by the jet grows 
from ^ 1 (at shock breakout) to ~ 3 in both models decreasing with latitude. 
At the end of the simulation Wmax is 29.35 (44.17) for model A (B), which is 
still smaller than the ones required for the fireball model ([||). However, our 
simulations have not been pushed far enough in time yet and, therfore, they can 
(at the present stage) neither account for the observational properties of GRBs 
nor of their afterglows. Instead, our set of numerical models can be regarded as 
simulations of a proto-GRB, because the scales treated in the simulations are still 
by more than 100 times smaller than the typical distances at which the fireball 
eventually becomes optically thin {^ 10^^ cm). 
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